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ITERATIVE GALERKIN DISCRETIZATIONS 
FOR STRONGLY MONOTONE PROBLEMS 

SCOTT CONGREVE AND THOMAS P. WIHLER 


Abstract. In this article we investigate the use of fixed point iterations to solve the Galerkin 
approximation of strictly monotone problems. As opposed to Newton’s method, which requires 
information from the previous iteration in order to linearise the iteration matrix (and thereby to 
recompute it) in each step, the alternative method used in this article exploits the monotonicity 
properties of the problem, and only needs the iteration matrix calculated once for all iterations 
of the fixed point method. We outline the abstract a priori and a posteriori analysis for the 
iteratively obtained solutions, and apply this to a finite element approximation of a second-order 
elliptic quasilinear boundary value problem. We show both theoretically, as well as numerically, 
how the number of iterations of the fixed point method can be restricted in dependence of the 
mesh size, or of the polynomial degree, to obtain optimal convergence. Using the a posteriori 
error analysis we also devise an adaptive algorithm for the generation of a sequence of Galerkin 
spaces (adaptively refined finite element meshes in the concrete example) to minimise the number 
of iterations on each space. 


1. Introduction 


In this paper we study Galerkin approximations of strictly monotone problems of the form: 

find uG X : A{u^v) =0 '^v G X. (1.1) 

Here, is a real Hilbert space, with inner product denoted by (-1 ■)x and induced norm ||a;|| = 
y {x,x)x- Furthermore, A : AxA—J-Risa possibly nonlinear form such that, for any u G X, 
the mapping v i—>■ A(u, v) is linear and bounded. Moreover, we suppose that A satisfies 
(PI) the strong monotonicity property 

A(u,u — v) — A(y,u — v) > co\\u — v\\‘j( \/u,v G X, (PI) 

for a constant cq > 0, and 
(P2) the Lipschitz continuity condition 

|A(u, w) — A(?;, w;)| < L||u — v||x||ii'||x 'iu,v,w G X, (P2) 

with a constant L > 0. 

Under these assumptions, there exists a unique solution u G X ot the weak formulation CH) ; see, 
e-g-; [HI Theorem 2.H] or [11]. In addition, the solution can be obtained as limit of a sequence 
it^,... G A resulting from the fixed point iteration 

{u'^,v)x = {u"~^,v)x - ^A{u'^~^,v) WvGX, n > 1, (1.2) 

for an arbitrary initial value G X. Indeed, defining the contraction constant 



there holds the a priori bound 

rn 

\\u-u-\\x<Y^\\u°-u^\\x, n>l, (1.4) 
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for the iteration (II.2L i.e., ||u — u^Wx 0. 

Restricting the iteration (ll.2p to a finite dimensional linear subspace Xh Q X, leads to an 
iterative Galerkin approximation scheme for (EH. More precisely, we consider, for an initial 
guess € Xh and n > 1, the iteration 

ulGXh-. {ul,v)x = {ul~^,v)x - yvGXh, (1.5) 

where cq and L are the constants from m and (IP2I) respectively. We emphasize that the problem 
of finding m)) from in the iteration scheme (EH) is linear and uniquely solvable. Similarly as 
for (11.21) and (II.IL the fixed point iteration (11.51) converges to the (unique) solution Uh S Xh of 
the Galerkin formulation 

A{uh,v) = {) Vv € Xh- (1.6) 

Furthermore, we note the a priori bound 

hn 

\\uh - ulWx < - uiWx, n > 1, (1.7) 

analogous to E3. 

In solving nonlinear differential equations numerically two approaches are commonly employed. 
Either the nonlinear problem under consideration is discretized by means of a suitable numerical 
scheme thereby resulting in a (finite-dimensional) nonlinear algebraic system, or the differential 
equation problem is approximated by a sequence of (locally) linearized problems which are dis¬ 
cretized subsequently. The latter approach is attractive from both a computational as well as an 
analytical view point; indeed, working with a sequence of linear problems allows the application 
of linear solvers as well as the use of a linear numerical analysis framework (e.g., in deriving error 
estimates). In the context of fixed point linearizations (11.51) yet another benefit comes into play: 
solving for from involves setting up and inverting a mass matrix on the left-hand side 

of ()1.5L We emphasize that this matrix is the same for all iterations; hence, it only needs to be 
computed once (on a given Galerkin space). 

The idea of approximating nonlinear problems within a linear Galerkin framework has been 
applied in a variety of works. For example in the article [S], the authors have considered general 
linearizations of strongly monotone operators, and have derived computable estimators for the total 
error (consisting of the linearization error and the Galerkin approximation error), with identifiable 
components for each of the error sources. A more specific linearization approach for monotone 
problems, which is based on the Newton method, has been presented in [7]. In a related context 
linear finite element approximations resulting from adaptive Newton linearization techniques as 
applied to seinilinear problems have been investigated in the papers m- Finally, we remark 
that the linear Galerkin approximation approach for nonlinear problems is not only employed for 
the purpose of obtaining linearized schemes, but also to address the issue of modelling errors in 
linearized models; see, e.g. mi- 

The aim of the current paper is to derive a priori and a posteriori error bounds for the Galerkin 
iteration method EH- Our error estimates are expressed as the summation of the linearization 
error resulting from the fixed point formulation with the Galerkin approximation error. In par¬ 
ticular, based on the a posteriori error analysis, we will develop an adaptive solution procedure 
for the numerical solution of EH that features an appropriate interplay between the fixed point 
iterations and possible Galerkin space enrichments (e.g., mesh refinements for finite elements); 
specifically, our scheme selects between these two options depending on whichever constitutes the 
dominant part of the total error. In this way, we aim to keep the number of fixed point iterations 
at a minimum in the sense that no unnecessary iterations are performed if they are not expected 
to contribute a substantial reduction of the error on the actual Galerkin space. 

The outline of the rest of this article is as follows. In Section [2] we derive an abstract analysis 
for the fixed point iteration EH, which includes the derivation of both a priori and a posteriori 
error bounds; in addition, we formulate an abstract adaptive procedure. The purpose of Section[3] 
is the application of our abstract theory to the finite element approximation of a second-order 
elliptic quasi-linear elliptic diffusion reaction boundary value problem; in particular, we derive a 
fully adaptive algorithm based on suitable a posteriori error estimates, and provide a series of 
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numerical experiments. Finally, in Section 2] we summarise the work presented and draw some 
conclusions. 


2. Abstract analysis 


2.1. Fixed point Galerkin approximation. As previously discussed, we let Xh be a finite 
dimensional linear subspace of a Hilbert space X. Then, in order to approximate the solution 
u € X oi ()l.ll) . we consider the Galerkin solution Uh S Xh defined in (ll.6|l . For the purpose 
of calculating Uh we consider, in turn, the discrete and linear fixed point iteration scheme (jl.5ll . 
Evidently, this is equivalent to a linear algebraic system of equations. More precisely, using basis 
functions (j)i £ A;^, for i = 1, ..., m, where m = dim(A/i) is the number of degrees of freedom, and 
letting 

m 

i=l 

for some unknown coefficients a" = obtain the linear system version of the fixed 

point iteration (11.51) : 


^ = E for J = 1,..., m. 

i—1 i—1 

Here, M is the iteration matrix defined by Mij = and with 


is the vector form of A(uh,v). We can see that the iteration matrix M does not depend on the 
iteration number n; hence, it only needs calculating once for all iterations of the fixed point method 
(on a given Galerkin space). 


2.2. A priori error bound. Denoting by 


el = u-ul (2.1) 

the error between the solution u of CH) and from dESl), we employ the triangle inequality 
and dni) to obtain 

un 

\\el\\x<\\u-un\\xAY^\\ul-ul\\x, 

where Uh £ Xh is the Galerkin solution defined in (11.61) . Furthermore, employing the monotonicity 
property lEU leads to 

co||m - Uh\\\ < A{u, u-Uh) - A{uh, u - Uh) 

= A(m, u — v) — A(uh,u — v), 
for any v £ Xh- Involving (IP2I) . we conclude 

Co||u - UhWx < ^11^ “ Uh\\x\\u - v\\x Vu £ Xh, 


and thus 

\\u-Uh\\x < —\\u-v\\x Vv £ Xh- 
Co 

Gombining these estimates we obtain the following result. 


Proposition 2.2. For the error between the solution u £ X of (HU and its iterative Galerkin 
approximation u]) £ Xh from (HU there holds the a priori error estimate 

Tj 

\\u-ul\\x<- inf \\u-v\\x + ^\\nl-ul\\x, 

Co v^Xh 1 — k 

for any n > 1. 
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2.3. A posteriori error analysis. In order to derive an a posteriori error analysis for (HU let 
us consider the auxiliary problem of finding u" G X such that 

Vu G X, n > 1. (2.3) 

We note that u” G X is a reconstruction (cf. [10]) in the sense that u)) G from (11.51) is the 
Galerkin approximation of u". We assume that we can bound the error between the solution u" G 
X of ()2.3I) and its Galerkin approximation G Xh in terms of an a posteriori computable 
quantity r]{ul,Xh), i.e., 

- KWx n>l. (2.4) 

Involving the monotonicity property (ED, the error e)) from (ED satisfies 

c,\\el\\\ < A{u,el) - A{ul,el) = -4l«,e)(). 


Furthermore, recalling (12.31) . we write 


c4el\\\<A{ul 




)-A{ 




) + —{u^-ul 
Co 


‘■h 




X 


= A{ul-\el) - Aiul^el) + — « - ul-\el)x + —(u" - <,e)()x. 


Then, using 


,n-l „n\ A(„,n „n\ , ^ ( n „,n-l n\ 

- \^h ~ I 

Co Co 

and applying the Gauchy-Schwarz inequality, we infer that 


Co 


colKrx < L\K-^ - KWxlKWx + —\\K - ul-^\\x\\el\\x + —||u" - K\\x\\el\\x, 


Co 


and dividing by co||e))||x, we obtain 


\e-h\\x < 


Co 


1 + 


Co 


\\K-<-^\\x + 




\U -u^Wx- 


Hence, inserting (12.41) . the following result can be deduced. 


Proposition 2.5. For the error between the solution u € X of dm and its iterative Galerkin 
approximation G Xh from dm there holds the a posteriori error estimate 

||e)(|U < - f 1 + -) \\< - K-^\\x + Xhl (2.6) 

Co V Co/ eg 

where r]{u'^,Xh) is given in (12.41) . 


2.4. An abstract adaptive algorithm. The a posteriori error estimate (12.61) shows that the 
error e)) from dm is controlled by two separate parts: a fixed point iteration error given by 
Lcq^ (1 + Lcq^) IK - u^-^Wx, and a Galerkin approximation term L^Cq Xh). When per¬ 

forming the fixed point iteration (11.51) it is worth noting that once the fixed point error is less than 
the Galerkin error carrying out another iteration will not cause a substantial reduction of the error 
on the actual Galerkin space. Based on this observation we are able to develop an algorithm which 
generates a sequence of hierarchically enriched Galerkin spaces Xhp C Xh ^2 C Xh^o C ... C X, 
with the aim of performing a minimal number of fixed point iterations at each enrichment step. 
Our algorithm will, therefore, feature an interplay between fixed point iterations and Galerkin 
space refinements. 

On the Galerkin space Xh^, * > 1, we define the Galerkin approximation error by 

^Galerkin.* •= —Vi^h^i, X h^i) , 

Cp 


and the fixed point error 





X- 


This allows us to write the a posteriori error bound (|2.6p as 


||u 


^Galerkin,2 


T £ 


n 

FP,2- 
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Here, we denote by uJJ ^ € Xfi,i the Galerkin solution obtained after n steps of the fixed point 
iteration (Hm on the current space X^y, for i > 0, the initial guess ^ G X^ i on the current 
Galerkin space Xh^i is obtained as the natural inclusion (or a projection) of the solution G 

Xh^i-i of the last (namely, the n*-th) iteration on the previous Galerkin space Xh,i-i to the space 
Xh^i- In particular, the fixed point iteration index n is reinitiated in each space enrichment step. 

Algorithm 2.7. Choose an initial starting space Xhp, and an initial guess q ^ Xih,o- 
for i •<— 0,1, 2,... do 
n ^ 0 
repeat 

n n + 1 

Perform a single fixed point iteration (US! to calculate itj) ^ G Xh,i. 
until fpp j < 'd^Galerkin,! 

Perform a hierarchical enrichment of X^ i 

based on the error indicator rjiu'^ Xh^i) from (12.41) 
to obtain a new Galerkin space X^ i^i D X^i. 

Define ^ ^ (by inclusion X^ i^i X^i or by projection) 

end for 

Here, d > 0 is a prescribed parameter. The algorithm is stopped if either the iteration number i 
reaches a given maximum, or if the right-hand side of (1^ is found to be sufficiently small. 

3. Application to quasilinear elliptic PDE 

3.1. Problem formulation. In this section, we focus on the numerical approximation of second- 
order elliptic diffusion reaction boundary value problems of the form 

—V • (/r(a;, |Vit|)Vit)-I-/(a;, u) = 0 in H, (3.1) 

M = 0 on r, (3.2) 

where H is a bounded, open, polygonal domain in with boundary F = dD. Here, we assume 
the following monotonicity conditions on the nonlinearities p, and /: 

(1) /i G C°(0 X [0, oo)), and there exist constants oi > 02 > 0 such that the following property 
is satisfied: 

a 2 {t — s) < p{x,t)t — p{x,s)s < ai{t — s), t>s>0, xGD. (3.3) 

(2) / G C°(H X R), and there exist constants fii > ^2 > 0 such that 

fi 2 it-s)< f[x,t)-f{x,s)<j5i{t-s), t>s, £c G H. (3.4) 

From [H Lemma 2.1] we note that, as p satisfies (|3.3I1 . for all vectors u, ta G R^ and all a; G O, we 
have 

\p{x, |u|)u — p{x, |ia|)ia| < q:i|u — ia|, (3.5) 

a 2 |u — < {pix, |u|)u — p{x, |ui|)ui) ■ (y — w). (3.6) 

Similarly, as / satisfies (13.4|) . it holds that for all s, t G R and all a; G H, 

\f{x,t)-f{x,s)\</3i\t-s\, (3.7) 

/32\t-s\^ < {f{x,t)- f{x,s)){t-s). (3.8) 

For ease of notation we shall suppress the dependence of p and f on x and write p(t) and f{u) 
instead of p{x,t) and f{x,u), respectively. 

The weak formulation of the boundary value problem (I3.ip - (I3.2I) is to find u G X := Hq{D) 
such that 

A(m,i;) = 0 yvGH({n), (3.9) 

where 

A{u,v) = j {p{\Xu\)Xu ■ VVf {u)v} dx, u,vGHq{D). (3.10) 

Jn 
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Throughout this section, we use function spaces based on a polygonal Lipschitz domain D C 
We denote by H^{D) the Sobolev space of order fc € Nq endowed with the norm || • In the 

case that fc = 0, we set H^{D) = L 2 {D) and denote the matching norm by || • Furthermore, 

we define Hq{D) as the space of functions in H^{D) with zero trace on dD. 

Introducing the inner product 

{u,v)n := / {a 2 ^u ■ Wv + 132 Uv} dx, u, G I/q(11), 

Jq 

where a 2 and /32 are the constants from (13.3|) and (13.4|) . respectively, we note the induced norm 

11^11^= «2||Vu||i^(^)+^2|kllL(n) 

on iJQ(n). 

3.2. Basic properties. Under the conditions (13.31) and (13.41) we can show that the properties (IPII) 
and (IP2I) are satisfied for X = i?g(U), and || • ||x := |||•|ll^• Indeed, noting the Poincare inequality, 

IkllL.(O) < C'pllV^^llL.(n) Vu G H^m, (3.11) 

where Cp > 0 is a constant dependent only on U, there holds the ensuing result. 


Proposition 3.12. Provided that (13.31) and (13.41) hold, then the form A from (13.101) is both strongly 
monotone with constant cq = 1 in (EU, and Lipschitz continuous with constant 


oi + max(/3i,“i/52/a2) Cl, 

- a2 + P2Cl 

in (IP2I) . Here, Cp > 0 is the Poincare constant from (13.111) . 

Proof. In order to show m with Co = 1, we apply (13.61) and (13.81) to arrive at 

A{u,u — v) — A{v,u — v) = / (^(|Vu|)Vu —/x(|Vu|)Vu) • V(u — u) da; 

Jn 

+ f {f{u) - f{v)) {u-v)dx 
Jn 


(3.13) 


> 


J |a 2 |V(u - v)f + I32\u - da; = |||u - uK 


Furthermore, to prove the Lipschitz continuity property 
combination with the Cauchy-Schwarz inequality, this yields 


we recall (ixa and (13.71) . In 


|^(u, w) — A(u, w)! < / |^(|Vu|)Vu —/i(|Vu|)Vu||Vw| da; + / \f{u) — f{v)\\w\dx 
Jn Jn 

< ai\\X{u - p)|lL 2 (n)ll^^llL 2 (n) + “ '^ 11 ^ 2 ( 0 ) 11 ^ 11 ^ 2 ( 0 )- 

We first consider the case when (32 = 0; then, noting that |||u|||n = we apply the 

Poincare inequality (13.111) to observe that 

\A{u,w) - A{v,w)\ < ai||V(u-u)||^^(f 2 )||Vu;|l^^(fj) +/SiCpIlV(u - || 

ai + /3iCp 


02 


||u-i^|||o|||w|||n. 


When /32 > 0 we introduce a constant 0 < 6 < (Ji and apply the Poincare inequality (13.111) , to get 
that 


\A{u,w) - A{v,w)\ < oi||V(m- + (/3i - (5 )||m- 

+ 5C|,||V(u-u)||^^(^)||Vu;|L2(0) 

+ ^ll^“^lli2(o)V^lkllL2(o)- 
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Using the Cauchy-Schwarz inequality yields 

\A{u,w) - A{v,w)\ < L((5)|||'u-i;|||n|||w|||n, 

where 

fai+SCl 

L{6) = max -, —- 

V “2 P2 

Minimizing L[5) within the given range, 0 < <5 < /3i, depends on the constants ai, 0 : 2 )/3i, /32- More 
precisely, if “i/as > then 5* = 0 is the optimal choice, and there holds L(0) = “i/as > 1; 
otherwise, we select 

/ 3 ia 2 — P2011 


S* ■= 


0-2 + 


S (0,/3i), 


and thereby obtain 

This completes the proof. 


L(<5*) = > 1. 

02 + (32Cp 


□ 


Remark 3.14. Incidentally, referring to, e.g., m Theorem 3.3.23] or [141 Theorem 2.H], the 
above result. Proposition 13.121 guarantees the existence of a unique solution u £ llo(U) of (13.91) . 

Remark 3.15. We note that the fixed point iteration 113 ) for the current problem (I3.1I) - (I3.2D 
reads in strong form as 

-a 2 AM" + = -a 2 An’"-i + /32u"-^ - (-V • + /(n"”^)) in O 

u" = 0 on dn, 

in the dual space of Ifg(U), for n > 1. 

3.3. Finite element discretization. In order to solve (13.91) by a fixed point Galerkin iteration, 
we will use a finite element framework. 


3.3.1. Meshes and spaces. We consider regular and shape-regular meshes Th that partition the 
domain U C into open disjoint triangles and/or parallelograms K, such that U = UifeTh 
We denote by Hk the elemental diameter oi K € Th, and let h = max^^gTy Hk- 

With this notation, for a fixed polynomial degree p > 1, we are ready to introduce the finite 
element space 

Ufem := {v G Hl{n) : v\k G Sp{K) VK G %} C (3.16) 


where 


Sp{K) 


Vp{K) if K is a triangle, 
Qp{K) if K is a parallelogram. 


Here, Vp{K) denotes the space of polynomials of total order at most p on K, while Qp{K) is the 
tensored space of polynomials of order at most p in each variable on K. 


3.3.2. Iterative Galerkin FEM. Based on the class of spaces Ufem introduced before, we can now 
introduce the finite element formulation for a linear fixed point formulation (13 of dUl): Given 
an initial guess G UfeM) we iterate for n = 1, 2,3..., 

{ul,Vh)Q = L~'^A{ul~^,Vh) yvheVfEM- (3.17) 

Remark 3.18. Recalling (II. 3p the contraction constant for this iteration is given by 

k=Vl- L-^ < 1- 

Here, we point out that, in the singularly perturbed case when 0:2 ~ ai = 0{e), for 0 < e <C 1, 
and (32 ~ Pi = 0(1), the contraction factor k does not deteriorate to 1. Indeed, this follows from 
the fact that the Lipschitz constant L from Proposition 13.121 remains robustly bounded from 0 in 
this situation. 


3.4. Error Analysis. We will now apply the abstract analysis derived in Section [2] to the iterative 
Galerkin method (13.171) for the numerical approximation of (|3.1I) - ()3.2D . 
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3.4.1. A priori error bound. Using our abstract a priori error analysis from Proposition 12.21 and 
applying suitable /ip-approximation results (see, e.g., [3]), we obtain a bound for the error between 
the numerical solution obtained at the n-th step of the fixed point iteration (13.171) and of the 
exact solution u from (1331) . For simplicity of presentation we assume a (quasi-) uniform diameter 
h > 0 of all elements. 


Theorem 3.19. Let u G H with k > 1, be the solution to the weak formula¬ 

tion dSH), G VpEM initial guess, and G Vfem the numerical solution after n steps of the 
fixed point iteration (13.171) ; then, for n > 1, there holds the a priori error estimate 

Lmin(K,p) , 

|||u - ulW^ < + 2L^ (1 - L-2) IK - rI K ’ (3-20) 

where C > 0 is a constant independent of h, p, and L from (13.131) . but depends on a 2 and /32 
from (1331) and (1331) . respectively. 

Remark 3.21. From the above Theorem 13.191 it is possible to predict the (approximate) number 
of fixed point iterations required to obtain an optimal convergence rate in the linear finite element 
iteration (13.171) . To this end, we ask for the second term on the right-hand side of (13.201) to converge 
at least at the rate of the first term. In order to discuss the resulting convergence behaviour of 
the numerical solution u)) obtained from (13.171) . we distinguish two different cases: 

• h-FEM: We fix a low polynomial degree p and investigate the convergence properties with 
respect to the mesh size h as h —>■ 0 (mesh refinement). Here, for k> p, we need 

(1 - L - 2)"/2 = O (hP) , 
and hence, n = O {\ log/i|) as h ^ 0. 

• p-FEM: We now fix the mesh, and suppose that the solution of (I3.1|) - (j3.2p is analytic. 
Then, as p —>■ oo, it can be shown that the FEM converges exponentially (see [12]), i.e., 
the error bound (|3.20p reads 

|||u - <|||„ < O (e-"^) + 2L^ (1 - UK - ^ 

for some constant b > 0. Again, balancing the two terms on the right, we require n = O {p) 
iterations as p —>■ oo. 

We will test these observations with some numerical experiments in Section 13.61 


3.4.2. A posteriori error analysis. In this section we obtain an a posteriori error bound for the 
error between the numerical solution uj) obtained at the n-th step of the fixed point iteration (13.171) 
and of the exact solution u obtained from (I3.1I) - (I3.2I) . According to our analysis in Section [2731 the 
key is to derive an a posteriori error estimate between the reconstruction u" G i7g(fl) from (12.31) 
and the iterative Galerkin solution from (USD (i.e., n)) from (13.171) in the present context); 
see (12.41) . 

To establish such a bound, we begin with a quasi-interpolation result. 

Lemma 3.22. Consider a finite element mesh Th, and a corresponding EEM space Vfem o,s 
in (13.161) . Moreover, let tt : Ed^iLl) Vfem be the Clement interpolation operator Then, 

E 

Ken 

for all V G £[^{^1), with a constant Cj > 0 independent of the local element sizes, and 

min(a^ih^,/3^i) if (32^0, 
otherwise, 



-V aaliV {v - -V “ kv\\ 


(OK) j ^ ^i\ 




for any K gTh- Here, a 2 and j32 are the constants from (13.31) and ([3.41) . respectively. 
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Proof. Let v £ Hq(Q,). We begin by recalling the following well-known approximation properties 
of the Clement interpolant: 


hj/\\v - TTvWl^^j^^ + l|V 

for any K G Th, with a constant C > 0 independent of the local element sizes and of v; for K G Th 
we denote by ujk the patch of all elements in Ph adjacent to K. In particular, following the 
approach in m. this implies that 

11 ^ - ^ + /32|klli,(,,^)) ^ 

and that, if /32 7 ^ 0 , then 

11^ - ’"^IIL(a) < («2|| -b /52||v|li,(,^^)) . 


Therefore, 

11^ - ^«llL(if) ^ («2||Vz;||^^(^^) -b /32|klli,(,^^)) > 

and so 

iK^h - + « 2 ||V (w - < C (a 2 ||Vz;||^^(^^) -b /32|k|li,(,^^)) • (3-24) 

Moreover, using the multiplicative trace inequality, that is, 

11 ^ - ^’^Wl^OK) < C' - ttv\\\^(^k) + ll^' - ^ 

we infer that 

h - ^^Wl^idK) < + Ik (a2||+/32|Mli,(^^)) ■ 

Observing that 


L-l I V 2 - 1/2 
n-K IK + 7 a «2 


_(1-1 1/2 1/2 
— [hx 0.2 Ik 


1/2 - 1/2 , „ 1 /; 

+ 1 j 7 a ^2 <^1k 


-/ 2 „,-V 2 


a 


we now arrive at 

a2'^'7A - ttiiIILoa) < C (a 2 || -b /32|k|li,(,^^)) • (3-25) 

Finally, combining p.24l) and (I3.25F and summation over all iF £ 7ii concludes the argument. 

□ 


In order to formulate the following result, we consider a series of meshes, {Th,i}i>o] for each 
mesh Th,i we denote the finite element space on that mesh as VpEM.i- 

Theorem 3.26. Let u £ be the exact solution to the boundary value problem (13.11) - (13.2F 

and Th,o be an initial mesh with initial guess q £ Vfem.o- Moreover, denote by Th,i the mesh 
after i mesh refinements, and let G Vfem.i be the FEM solution obtained after n steps of the 
fixed point iteration (j3.17|) on Th,i. Here, the initial guess ^ £ VpEM.i on the current mesh Th,i, 
* > 0, is obtained as an (appropriate) projection of the solution ujj G VpEM.i-i of the last (n*- 
th) iteration on the mesh Th,i-i to the space Vfem,!- Then, for n > 1, there holds the a posteriori 
error estimate 

/ \ V 2 




E 


+ i (1 "b L) 


. KeTh,: 


‘‘h,i 


n-1 

^h,i 


where Cj is the constant from Lemma \3.221 and 


Vk = PK 


( 


Vu 


n— 1 
h.i 


-V 2 _^^V 2 


Vm 


)v<-^)+L^F 
)Vul-^+L^a2V (<,-<-') 


L2iK) 


L2{dK\r) 
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for any K G Th,i o-nd n> 1, are local error indicators. Here, ■jk, K G Th,i, is defined in (|3.23l) 
and 


F(p) = —a2^v + 132V. 

Moreover, for an edge e C dK'^ fl dK~ between two neighbouring elements G Th,i, we signify 
by lnJlg = v'^\e • viii+ + v~\e ■ the jump of a (vector-valued) function v along e, where 

denote the traces of the function v on the edge e taken from the interior of K^, respectively, and 
nx± are the unit outward normal vectors on dK^, respectively. 

Proof. Recalling our abstract result, ProDOsition l2.5l it is sufficient to derive a quantity r]{u]) Vfem.i) 
such that 


II’W" - K.ilin < ■n{uh,i, Vfem.O; 


cf. (12.41) . The reconstruction it" G fulfills 


«,i;)n = {ul/,v)n-L ‘^Aiv 


'iv G D Vfem.i, 


see (12.3L where u)) ^ G Vfem.i is just the Galerkin approximation of it" (cp. (13.171) 1. 

Define the error ejj j = u" — itjjj, and let Vh,i = i ^ Vfem.z, where tt is the interpolation 
operator from Lemma [3.221 We notice that there holds 


lll~n |||2 \ f^.Ti 'ffn \ 

= - V,,,) - (<, - . 


Integration by parts elementwise leads to 




KeTh 


KdTh.i'’^ 

E (/ {<f) + (<, - <-')) - VhA Ax 

) (ul^ - ulf)) {Ai, - vhA dx 

E (/ (<, - KA)) - i;,,) da: 

a2V (ul^i - uIA)) ■ nK {Ai,i - Vh,i) ds. 


K&Th.i 


KeTh 

- ^ LM 


n—1 

h,i 




A few elementary calculations show that 


s LM 


Vu 


n—1 

h.i 


Vulf + L^a2V (ul, - ulf)) ■ riK ^Ki - Vh^ ds 




K^Th,, 
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and thus, using the Cauchy-Schwarz inequality, implies 


ATeTh,, 


L2{K) 


^lK"\K^-yhA\L,^,K) 


: E 

K&Th,i 


-1/4 1/4 
^2 Ik 






1/2 


n ^ 

\K^Th,i 


Vk 


V 2 


.K(iTh,i \ 


L 2 (a-R\r) 


iL 2 (aAr\r) 


1/2 


Therefore, we infer that, employing Lemma 13.221 


which implies 


V ifGTfe.i / 

J,ill|n < C/L”^ f ^ 77!^^ =: Vfem.j)- 

V AGTfc.i / 


Inserting this bound into (12.611 with cq = 1 (cp. Proposition 13.1211 completes the proof. □ 


3.5. Adaptive refinement algorithm. Proceeding along the lines of Sectionwe notice that 
the a posteriori error bound from Theorem 13.261 controls the error in terms of two contributions: 
The finite element error, defined as 


£ 


n _ 

FEM,2 — 



1/2 


and the fixed point error 


^FP,z — -^(1 + 



— u 


n— 

h,i 


1 


This allows us to write the error bound as 


a' 


\U — Ui, 


<Ci£ 


FEM, 


-FP,; 


Based on this bound we can cast the abstract adaptive Algorithm 12 . 71 into the fixed point Galerkin 
iteration (13.1711 for the solution of (I3.1I1 - (I3.2I1 . 


Algorithm 3.27. Choose an initial starting mesh Xifl, and an initial guess u® g € Vfem.o the 
associated finite element space Pfem.o (of fixed polynomial degree p > 1). 

for i ^ 0,1, 2 ,... do 

n ^ 0 
repeat 

n n + 1 

Perform a single fixed point iteration (13.1711 to calculate ^. 

until 

until Cpp^ ^ '^^FEM.i 

Perform mesh refinement (and/or derefinement) on Th,i 

based on the error indicators rjK from Theorem \3.26\ 

together with a suitable marking strategy in order to obtain Xh^i+i- 

<.i+l ^ 

end for 

Here, TTi^i+i is some projection from Pfem.i to VFEM.i+i (for instance, the /^-projection), and 
i9 > 0 is a (prescribed) parameter. 
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(a) (b) 


Figure 1. Convergence with restricted nnmber of fixed point iterations with nniform (a] 
p-refinement and (b) /i-refinement. 


3.6. Numerical experiments. In this section we perform a series of numerical experiments 
to validate the a priori and a posteriori error bounds for the fixed point iteration (j3.17l) from 
Theorem 13.191 and Theorem 13.261 respectively. 


3.6.1. Validation of E,emark \S.21\ We consider the domain Fl = (0,1)^ C with nonlinearity 

m(|Vu|) = 2+ ^ 

1 + I Vm| 

and select / independent of u such that the analytical solution to (I3.1I) - (I3.2I) is given by 

u[x,y) = x{l - x)y{l - y){l - 


We note that = /?2 = 0, oi = 3 and 02 = 

Firstly, we consider the case when the mesh is fixed as a 16 x 16 uniform square mesh of quadri¬ 
laterals and perform uniform refinement of the polynomial degree p from an initial guess q = 0 . 
In this situation we restrict the number of iterations of the fixed point iteration to Cp ■ p, for 
Cp = 1,2,3 and plot in Figure [ Jl^a)| the error |||u — u/t|||n against the polynomial degree p. For 
comparison, we also perform the same experiment continuing the fixed point iteration until the 
residual A{u]^,v) is below a given tolerance (10“^"^) and, hence, the approximation is close to the 
best possible FEM approximation for the mesh. We clearly observe that by restricting the number 
of iterations we obtain exponential convergence of the error, and when Cp = 3 we gain the same 
convergence rate as allowing the iteration to continue until a tolerance is reached. Hence only 
performing the iteration 3p times gives an optimal convergence rate in the given example. 

Secondly, we consider a fixed polynomial degree p = 2 and perform h refinement to generate a 
sequence of 2^ x 2^ uniform square meshes of quadrilaterals, for iV = 3,..., 8 . We again perform 
both a restriction of the number of iterations of the fixed point method to Cn ■ N, for Cn = 1 , 2 , 
as well as allowing the iteration to continue until a tolerance is reached, and plot in Figure [ ^b)| 
the error |||u — u/j|||n against | log 2 h\. We obtain algebraic convergence, and already when Cn = 2 
we achieve the optimal convergence rate 0{h^). 


3.6.2. Validation of Theorem \S.2fi\ and E,emark \S.18\ We now consider automatic h-adaptive mesh 
refinement, with linear (p = 1 ) basis functions, using Algorithm 13.271 and the a posteriori error 
bound from Theorem 13.261 For the purpose of mesh refinement we use a fixed fraction refine¬ 
ment strategy, where the 25% of elements with the largest local error indicators pK are marked 
for refinement, and the 5% of elements with the smallest local error indicators are marked for 
derefinement. 
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(c) 



(d) 


Figure 2. Example 1. (a) Error in |||•|||n-norm and error bound from Theorem 13.261 


after the final fixed point iteration on each mesh compared to the number of degrees of 
freedom; (b) Effectivity at each fixed point iteration for all meshes; (c) Number of fixed 
point iterations on each mesh; (d) Mesh Th,! after 7 /i-refinements. 


Example 1. Nonlinear diffusion We first consider the case of a u-independent / with a nonlinear 
pL on the unit square = (0,1)^ C To this end, we let 

fi{\yu\) = 1 + arctan(|VM|^), 

and select / such that the analytical solution to (EH)-(1321) is given by 

u{x, y) = x{l - x)y(l - y)(l - 

We note that /3i = /?2 = 0, cti = 1 + v ^/2 + and a 2 = 1- For this problem we set the steering 
parameter d = 1/2 in Algorithm 13.271 

We first plot, in Figure Qfa) the relative true error error bound, with 

Cl = 1, from Theorem l3.26l after the last iteration n* on each mesh i against the number of degrees 
of freedom on that mesh. As can be seen, both the true error and the error bound converge at 
the same rate, and the error bound appears to overestimate the true error by a roughly constant 
amount. We also consider in Figure [ ^b)| the effectivity index at each step of the fixed point 
iteration on each mesh, where the effectivity index is the error bound (calculated with C/ = 1) 
divided by the true error. As can be seen this is roughly constant (approximately 4) for all meshes 
and iterations, which indicates that the error bound overestimates the true error by roughly this 
amount, independent of mesh properties. We do note, however, that the effectivity rises slightly 
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(c) 



Figure 3. Example 2. (a) Error in |||•|||n-norm and error bound from Theorem 13.261 


after the final fixed point iteration on each mesh compared to the number of degrees of 
freedom; (b) Effectivity at each fixed point iteration for all meshes; (c) Number of fixed 
point iterations on each mesh; (d) Mesh Th ,7 after 7 /i-refinements. 


due to the fixed point iteration on each mesh, this is likely caused by setting C/ = 1. We also 
plot in Figure l ^c) I the number of fixed point iterations at each mesh step required to ensure that 
the fixed point error is less than the finite element error. We note this is fairly constant although 
minor variations exist. A few mesh refinements are made early on which is likely caused by the 
fact that the features of the solution are not accurately captured at the beginning. 

In Figure l ^d) I we plot the mesh 7^ 7 after 7 h-adaptive mesh refinements. The areas of mesh 
refinement coincide with the hill and valley in the analytical solution, which is the location we 
would expect the greatest error to occur, and matches the sort of refinement that occurs when the 
nonlinear methods are computed to a minimal residual. This suggests that the mesh refinement 
algorithm behaves in the expected manner. 

Example 2. Strong nonlinear reaetion We now consider a fairly strong nonlinear / with a 
constant diffusion coefficient /r(|Vu|) = e, where e is a small positive constant, on the unit square 

= (0,1)^ C To this end, we consider e = 0.01 and let 


f{u) = i0.2 + x^ + y^) 


+ 1 


Fc[x,y), 
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where c(a:, y) is a function dependent only on x and y selected such that the analytical solution to 
(IXTD (15^ is given by 

2 /) = (1 - x){l - - l)(e''^^ - 1 ). (3.28) 

We note that ax = a 2 = £, /3i = and /32 = V®- For this problem we set the steering 

parameter d = 1 in Algorithm 13.271 

We again plot, in Figure the relative true error |||“-“h.i|||jj/|||tt|||t 2 and error bound, with 

C/ = 1, from Theorem 13.261 after the last iteration n* on each mesh i against the number of 
degrees of freedom on that mesh. Except for a few early meshes, where the mesh is unlikely to 
accurately capture the boundary layer in the analytical solution, both the true error and the error 
bound converge at the similar rate with the error bound overestimating the true error by a roughly 
constant amount. This is supported by the effectivity indices at each iteration, Figure [ ^fb)] which 
are roughly constant for all meshes and iterations, although slightly decreasing over the course of 
the mesh refinement. For this problem we note that the number of fixed point iterations at each 

mesh step is fairly high, caused by the stronger nonlinearity, but only once the boundary layer is 

captured accurately by the mesh (the 5-th mesh onwards); cf. Figure in Before this mesh the 
finite element error is considerably larger than the fixed point error due to the inaccurate capture 
of the boundary layer. The mesh after 7 h-adaptive mesh refinements. Figure [ ^^d)] demonstrates 
how the mesh captures the boundary layer. This demonstrates an important benefit of only 
refining the fixed point error while it is greater than the finite element error, as the algorithm 
has managed to reduce the number of iterations in the early meshes by a considerable number by 
not performing fixed point iterations until after the mesh has started to accurately capture the 
solution’s features. 


Example 3. Nonlinear reaction We now consider a weaker nonlinear / with a constant diffusion 
coefficient /i(|VM|) = e, where e is a small positive constant, such that /3i ~ /32 = 0{1), on the 
unit square = (0,1)^ C To this end, we consider e = 10“^, for fc = 0,..., 6 , and let 


f{u) = 


c{x,y), 


lOu'^ + l 

where c(a:, y) is a function dependent only on x and y selected such that the analytical solution to 
mi es is given by (I3.28p . We note that ai = a 2 = e, /3i = ss/so, and P 2 = 1- For this problem 
we set the steering parameter -d = 1 in Algorithm 13.271 

We plot, in Figure the true error u — for e = 10“^, where /c = 0,..., 6, on each 

mesh i against the number of degrees of freedom on that mesh. We note that we appear to achieve 
a higher initial rate of convergence for smaller e values, although they all appear to tend to similar 
convergence rates as refinement progresses. For each e we also calculate the effectivity indices 
at each iteration, Figure gKb)| We note that initially they are highly oscillatory for small values 
of e, but as refinements progress they tend to smoothen towards a constant value, with the high 
values of e converging at earlier mesh steps; in particular, the effectivity indices do not deteriorate 
as £ —)► 0+. We also plot in Figure in the number of fixed point iterations at each mesh step 
required to ensure that the fixed point error is less than the finite element error. We note this is 
fairly constant and independent of the value of £, which supports Remark 13.181 


4. Conclusion 

In this article we have shown that it is possible within a finite element framework to use a 
simple fixed point iteration to solve strongly monotone quasi-linear elliptic PDEs requiring only 
the computation of the iteration matrix, opposed to a Newton’s method requiring computation 
at each iteration. We have shown that an optimal a priori convergence rate can be obtained 
for a fixed number of iterations, dependent on the mesh size or polynomial degree. We have 
demonstrated that it is possible to perform adaptive mesh refinement based on an a posteriori 
error analysis, where only a minimal number of fixed point iterations are required on each mesh 
to obtain a good approximation to the solution, without continuing the fixed point iteration until 
the fixed point error is insignificant. 
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0 5 10 15 
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Figure 4- Example 3. (a) Error in |||•|||^-norm after the final fixed point iteration on each 
mesh compared to the number of degrees of freedom for various values of e; (b) Effectivity 
at each fixed point iteration for all meshes; (c) Number of fixed point iterations on each 
mesh. 
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